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ABSTRACT 

The launching of magnetohydrodynamic outflows from accretion disks is considered. We formulate 
a model for the local vertical structure of a thin disk threaded by a poloidal magnetic field of dipolar 
symmetry. The model consists of an optically thick disk matched to an isothermal atmosphere. The disk 
is supposed to be turbulent and possesses an effective viscosity and an effective magnetic diffusivity. In 
the atmosphere, if the magnetic field lines are inclined sufficiently to the vertical, a magnetocentrifugal 
outflow is driven and passes through a slow magnetosonic point close to the surface. We determine how 
the rate of mass loss varies with the strength and inclination of the magnetic field. In particular, we 
find that for disks in which the mean poloidal field is sufficiently strong to stabilize the disk against 
the magnetorotational instability, the mass loss rate decreases extremely rapidly with increasing field 
strength, and is maximal at an inclination angle of 40-50°. For turbulent disks with weaker mean fields, 
the mass loss rate increases monotonically with increasing strength and inclination of the field, but the 
solution branch terminates before achieving excessive mass loss rates. Our results suggest that efficient 
jet launching occurs for a limited range of field strengths, and a limited range of inclination angles 
in excess of 30°. In addition, we determine the direction and rate of radial migration of the poloidal 
magnetic flux, and discuss whether configurations suitable for jet launching can be maintained against 
dissipation. 

Subject headings: accretion, accretion disks — galaxies: jets — hydrodynamics — ISM: jets and 
outflows — magnetic fields — MHD 



1. INTRODUCTION 

Jets and other outflows are commonly observed from 
young stellar objects, interacting binary stars and active 
galactic nuclei. It is generally believed that the outflows 
are produced by the accretion disks in these systems. The 
widely differing properties of the central objects in these 
systems suggest that a mechanism may be at work that 
is largely independent of the nature of the central object 
(e.g. Livio 1997). 

A particularly promising mechanism for the acceleration 
and collimation of outflows involves a large-scale poloidal 
magnetic field that threads the disk. The influential model 
of Blandford & Payne (1982) established the significance 
of the angle of inclination, i, of the poloidal magnetic field 
lines to the vertical at the surface of the disk. When 
i > 30°, matter that is just above the surface, being forced 
to corotate with the foot-point of the field line, is acceler- 
ated outwards along it by the centrifugal force. At greater 
distances, the flow may collimate into a jet by magnetic 
hoop stresses or poloidal collimation. A clear review of the 
physics of magnetocentrifugally driven outflows has been 
given by Spruit (1996). 

Several groups have performed axisymmetric numerical 
simulations of outflows using this mechanism (e.g. Ustyu- 
gova et al. 1995; Romanova et al. 1997; Ustyugova et al. 
1999; Ouyed & Pudritz 1997a, 1997b, 1999; Krasnopol- 
sky, Li, & Blandford 1999). Like the model of Blandford 



& Payne (1982), these calculations do not resolve the verti- 
cal structure of the disk but treat it as a boundary surface 
that loads mass at a specified rate on to the magnetic field 
lines. While these simulations have convincingly demon- 
strated some aspects of the magnetocentrifugal model of 
jet production, several fundamental issues remain to be 
resolved. In this paper, we will attempt to address two of 
these questions. First, what determines the rate of mass 
loss in the outflow? How does this depend on the strength 
and inclination of the magnetic field, or on other proper- 
ties of the disk? Second, what is the long-term evolution 
of a large-scale magnetic field in a disk? Is it possible to 
assemble a magnetic configuration suitable for jet launch- 
ing, or to maintain it against dissipation? A proper under- 
standing of these issues is essential if we are to explain the 
conditions that regulate the production of astrophysical 
jets. 

The magnetocentrifugal mechanism can explain the ac- 
celeration of outflows even when the matter is 'cold' in 
the sense that the temperature is much less than the virial 
temperature or escape temperature. However, as men- 
tioned by Blandford & Payne (1982) and calculated in 
detail by Ogilvie (1997) and by Ogilvie & Livio (1998; 
hereafter Paper I), some amount of thermal assistance is 
still required to overcome the potential barrier between the 
surface and the slow magnetosonic point ('sonic point'). 
The barrier occurs because the field lines are not straight 
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inside the disk and the angular velocity deviates slightly 
from the Keplerian value because of the radial Lorentz 
force associated with the bending of the field lines (Shu 
1991; Wardle & Konigl 1993). A proper calculation of 
the vertical structure of the disk, including these effects, 
is therefore required in order to determine the height of 
the potential barrier and the rate of mass loss in the out- 
flow. Indeed, it is obvious from the well-known properties 
of transonic outflows that the mass loss rate depends on 
the physical conditions below the sonic point, and that the 
vertical structure must therefore be resolved. 

In Paper I we carried out such a calculation for disks 
that are rather strongly magnetized in the sense that the 
magnetorotational instability, which leads to turbulence in 
accretion disks (Balbus & Hawley 1998), is suppressed or 
nearly so. We assumed that the magnetic field enforces 
strict isorotation and showed that the potential barrier 
increases very steeply as the field is made stronger and 
the disk becomes more sub-Keplerian. This effect would 
suppress outflows from strongly magnetized disks unless 
an additional source of energy, such as coronal heating, 
were present (in accord with the suggestion made by Livio 
1997). 

A disadvantage of the calculation in Paper I is that no 
explanation was given for the source of the effective vis- 
cosity of a disk in which the magnetorotational instabil- 
ity is suppressed. The solutions did not extend to field 
strengths much below the stability boundary. When the 
field is weaker than this, it cannot be expected to enforce 
isorotation, and the model requires some modification. In 
addition, the effect of the turbulence on the mean mag- 
netic field should be modeled, most simply through the 
introduction of an effective magnetic diffusivity. 

However, the presence of a turbulent diffusivity may 
cause problems for the magnetocentrifugal mechanism. 
The analysis by Lubow, Papaloizou, & Pringle (1994a) 
suggests that, if the effective magnetic Prandtl number of 
the disk is of order unity, as might be expected, it is impos- 
sible to sustain a steady configuration with a significantly 
inclined field in a thin disk (see also Heyvaerts, Priest, 
& Bardou 1996; Reyes-Ruiz & Stepinski 1996). Although 
the accretion flow tends to drag magnetic flux inwards, 
the turbulent diffusivity expels flux faster if the inclina- 
tion is large. The inclination angle in a steady state is then 
expected to be comparable to the angular semi-thickness 
H/r of the disk. This constraint can be avoided if the 
accretion flow is due primarily to the loss of angular mo- 
mentum in an outflow, and might also be relieved by a 
dynamo operating in the disk if special conditions are met 
(Campbell, Papaloizou, & Agapitou 1998). 

The purpose of the present paper is to explore a model 
for the vertical structure of magnetized disks. We will ex- 
tend the analysis of Paper I to allow for the possibility 
of more weakly magnetized disks in which there is turbu- 
lence, and strict isorotation does not hold. We will also 
give a better treatment of the matching between the disk 
and the atmosphere by applying photospheric boundary 
conditions, allowing us to calculate the mass loss rate in 
the outflow explicitly. At the same time, we will solve 
for the rate of dragging of magnetic flux and thereby re- 
fine the previous estimates which may have been based on 
oversimple arguments. 

Related calculations of the vertical structure of mag- 



netized disks, and of the disk-jet connection, have been 
made by Konigl (1989), Wardle & Konigl (1993), Li (1995, 
1996), Ferreira (1997), Campbell (1999), Casse & Ferreira 
(2000), and Shalybkov & Riidiger (2000). Our method and 
results are significantly different from all previous calcula- 
tions and comparisons of the relevant issues will be made 
towards the end of this paper. 

2. CONTEXT OF THE PRESENT CALCULATION 

To set this calculation in its proper context, we con- 
sider here some of the wider issues relating to magnetized 
accretion disks and outflows. 

For a thin disk without a large-scale magnetic field, 
there is a clear division of the physical problem into two 
aspects. First, a model is required for the local vertical 
structure at any given radius r and time t. Such a model 
takes the surface mass density E as a parameter and pre- 
dicts quantities such as the vertically integrated viscous 
stress Q. The vertical structure may be assumed to be 
stationary on the dynamical time-scale. Second, the one- 
dimensional conservation equations for mass and angular 
momentum determine how E(r, t) evolves on the viscous 
time-scale, given the local relation between Q and E (e.g. 
Lynden-Bell & Pringle 1974). 

For a disk with a large-scale poloidal magnetic field but 
no outflows (e.g. a polytrope with a force-free or vacuum 
exterior), the situation is similar but more complicated. 
Now the local vertical structure at a given radius depends 
not only on E but also on the vertical magnetic field B z 
and the inclination angle i. In addition to an evolution- 
ary equation for E, there is a one-dimensional conservation 
equation for magnetic flux, which determines how the flux 
function ?/>(r, t) evolves on the viscous time-scale. B z is 
simply related to the radial derivative of ip. Finally, the 
force-free magnetic equilibrium in the exterior of the disk 
leads to a global relation between ip and the inclination 
angle i. The full problem therefore involves an integro- 
diffcrcntial equation. This has been considered by Lubow 
et al. (1994a), although they did not examine the effect of 
the magnetic field on the vertical disk structure. 

When outflows occur along the magnetic field lines, fur- 
ther couplings exist. The exterior magnetic field is modi- 
fied from the force-free solution by the inertial forces asso- 
ciated with the outflow. This effect depends on the amount 
of mass loading, and is therefore coupled to the launching 
problem. The exterior field is no longer purely poloidal but 
becomes significantly twisted beyond the Alfven surface, 
although non-axisymmetric (e.g. kink) instabilities may 
set in here (e.g. Spruit, Foglizzo, & Stehle 1997). This 
part of the solution is established on the Alfven travel 
time to the Alfven surface, which is comparable to the dy- 
namical time-scale of the disk. Finally, the mass loss and 
especially the angular momentum loss in the outflow feed 
back into the evolutionary equations for the disk. 

The full problem therefore involves three aspects: the 
local vertical structure of the disk and the launching prob- 
lem; the global structure of the exterior magnetic field and 
the dynamics of the outflow along those field lines; and 
the evolution of mass, angular momentum, and magnetic 
flux on a longer time-scale. The first aspect is the subject 
of this paper, while the second aspect is well described 
by some of the numerical simulations referred to above. 
The third aspect has been considered implicitly in steady 
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in dragging in magnetic flux. This can be explained by 
noting that the effective magnetic Reynolds number of the 
accretion flow, based on the disk thickness H, is small, of 
order H/r (Heyvaerts et al. 1996). Similar arguments sug- 
gest that a configuration in which the field lines are bent 
significantly from the vertical (e.g. to achieve i > 30°) 
cannot be sustained on the viscous time-scale. 

These arguments are based on a kinematic analysis of 
the magnetic induction equation in which the radial ve- 
locity and magnetic diffusivity are prescribed quantities. 
They also depend on a simple order-of-magnitude treat- 
ment of the vertical structure of the magnetic field. In 
this paper we will present a numerical treatment of a set 
of equations, including the induction equation, in which 
the radial velocity is self-consistently determined. This is 
important because, in a jet-launching configuration, the 
magnetic field must become dynamically dominant above 
some height and will then control the radial velocity. We 
will find results that differ significantly from the estimates 
of Lubow et al. (1994) under some circumstances. Before 
presenting the numerical model, we therefore reconsider 
the problem of the induction equation from an analytical 
viewpoint. 

We assume that the mean poloidal magnetic field is 
axisymmetric and may be described by a flux function 
ip(r,z, t) such that 



B r 



I dip _ 1 dtp 

r dz' z r dr 



(1) 



models such as that of Casse & Ferreira (2000), but the 
development of a complete evolutionary scheme for mass, 
angular momentum and magnetic flux in the general, non- 
steady case remains a challenge (see Lovelace, Newman, & 
Romanova 1997 for a simplified version of such a scheme). 4 

There has been some confusion in the literature regard- 
ing the interplay between these aspects of the problem. For 
example, studies of the launching problem have often im- 
posed the condition that the poloidal magnetic flux should 
not migrate radially (e.g. Konigl 1989; Li 1995). That is, 
there should be an instantaneous balance between inward 
advection of flux by the accretion flow, and outward trans- 
port due to turbulent, Ohmic, or ambipolar diffusion. This 
may not always be appropriate because the migration oc- 
curs on the viscous time-scale, whereas the solution of the 
launching problem need only be stationary on the dynam- 
ical time-scale. In Wardle & Konigl (1993) the rate of flux 
migration was treated as a free parameter. 

Moreover, some of these studies appear to predict the 
value of the toroidal magnetic field at the surface of the 
disk, B^ s , even when the dynamics of the outflow beyond 
the sonic point has not been considered (e.g. Campbell 
1999) or when there is no outflow (e.g. Shalybkov & 
Riidiger 2000). This is unsatisfactory because B^ deter- 
mines the external magnetic torque acting on the disk, 
which certainly depends on the dynamics of the outflow 
in the trans- Alfvenic region. These incorrectly specified 
studies may be motivated by the fact that the model of 
Blandford & Payne (1982) appears to require B^ to be 
determined by the disk (through their parameter A). How- 
ever, this is misleading and may be attributed to the fact 
that the Blandford & Payne model is underconstrained be- 
cause it is missing a boundary condition at large distances 
from the disk (cf. Ostriker 1997). Indeed, solutions of their 
model generically behave unphysically at large distances 
without having passed through the modified fast magne- 
tosonic point. In contrast, Krasnopolsky et al. (1999) have 
given a clear discussion of which quantities may or may not 
be specified as boundary conditions when treating the disk 
as a boundary surface. In their model B^, s is determined 
by the outflow, not by the disk. 

In our local study, therefore, we will not impose the 
constraint of zero flux migration; rather, we will impose 
the value of B,p s and determine the rate of flux migra- 
tion. Usually we will specify B<p s = 0, meaning that 
the outflow is absent or exerts only a weak torque. Ef- 
ficient magnetocentrifugal outflows are expected to have 
I -60s | "C \B Z \, with the toroidal field becoming comparable 
to the poloidal component only at greater distances from 
the source, where the Alfvcn surface is located (Spruit 
1996). 

3. THE EVOLUTION OF MAGNETIC FLUX 

As noted above, an important investigation of the evo- 
lution of the poloidal magnetic flux was made by Lubow 
et al. (1994a), who concluded that, if the disk is turbu- 
lent, with an effective magnetic Prandtl number of order 
unity, the accretion flow will be almost entirely ineffective 

4 There have been many calculations simulating the time-dependent development of a jet from a numerically resolved disc. For example, 
Kudoh, Matsumoto, & Shibata (1998) set up a thick torus with a non-rotating corona, and introduced a vertical magnetic field. Owing to the 
lack of equilibrium in the initial condition, a rapid adjustment occurs, and Kudoh et al. followed the evolution for only a single orbit or so. 
Unfortunately there is little reason to expect the properties of such transient phenomena to agree with quasi-steady models of jets from thin 
discs. 



For a thin disk containing a bending poloidal field of dipo- 
lar symmetry, the flux function has the form (Ogilvie 1997) 



ijj = ip (r, t) +^ 1 (r,z,i), 



(2) 



where \tpi\ -C IV'ol- Then B r and B z are comparable in 
magnitude if |V'i/V'o| — 0(H/r). The flux content of the 
disk is determined essentially by tpo, while -01 allows for 
bending of the field lines within the disk. 

We assume that the mean magnetic field in a turbulent 
disk may be treated within the framework of mean-field 
electrodynamics (e.g. Moffatt 1978). In the absence of a 
mean-field dynamo effect (or in the absence of a toroidal 
field), the flux function then satisfies the mean induction 
equation, 



~dt 



u-Vip — rjr V- 



(3) 



where u is the mean velocity and r] the turbulent magnetic 
diffusivity. In a thin disk the dominant terms are 



dt 



dr 



d 



^ dr \r 



dr 



+ V 



d^h 
dz 2 



(4) 



The neglected terms are all smaller than those retained 
because \ipi \ <C IV'o and \u r \ 3> \u z \. The first term on the 
right-hand side is also sometimes neglected (e.g. Lubow 
et al. 1994a). 
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Although this equation appears to be an evolutionary 
equation for ip it should be observed that the equation 
is defined for all values of z whereas ipo is a function of r 
and t only. Moreover, the equation is not closed because 
of the term involving tpi. Following Lubow et al. (1994a), 
we might attempt to close the equation by vertical aver- 
aging. Conventionally in accretion-disk theory one uses 
density-weighted averages u r and f\ defined by 



pu r dz, £77 



pndz, 



where p is the density and 

£ = J pdz 



(5) 



(6) 



is the surface density, the integrals being over the full ver- 
tical extent of the disk. The density-weighted average u r 
is particularly appropriate because it is directly related 
to the mass accretion rate. The vertical average of the 
induction equation is then 



chpo 
dt 



dr 



dr \r dr 



PH—dz. (7) 



Lubow et al. proceeded by approximating the last term as 
— (r/H)fjB rs . This closes the equation because B rs , which 
is proportional to the vertically integrated toroidal elec- 
tric current, can be globally related to the flux function 
by considering the force balance exterior to the disk. If 
the exterior region is treated either as an insulating vac- 
uum, or as a force-free medium with vanishing B^, the 
exterior poloidal field is potential. This leads to a global 
relation of the form 



B r 



(8) 



where £ is a certain linear integral operator (see also 
Ogilvic 1997). 

This order-of-magnitude estimate of the final term ap- 
pears reasonable if the field lines bend mainly in the dens- 
est layers of the disc near the midplanc. However, the term 
might be significantly overestimated if the bending occurs 
mainly in the upper layers of the disk where pi] is likely 
to be considerably smaller. Another way of averaging the 
induction equation reinforces this concern. If the equation 
is divided by r\ and integrated with respect to z, we obtain 



where 



dipo , dip d (I dtp \ 

at or or \r or ' 



'■*(r,t) = Jjdzf J^-dz, 
V*(r,t) = J dz I j ^dz, 



(9) 

(10) 
(11) 



and 



Q{r,t) = 2rj J ± dz. (12) 

Now the precise shape of the field line has truly been elim- 
inated in favor of B rs . This therefore appears to be the 



'correct' way of averaging the equation. However, we are 
now faced with vertical averages of u r and 77 weighted by 
I/77, not by p. Even if 77 is independent of z, these could be 
significantly different from the density-weighted averages, 
leading to different conclusions about flux dragging. 

Furthermore, as mentioned above, the presence of the 
magnetic field will in general change the vertical profiles 
of u r and n. If the field is extremely weak so that the in- 
duction equation may be treated 'kinematically', then the 
1/77 averaging method is surely the correct one. However, 
the vertical profile of u r , in particular, is determined by 
subtle effects (e.g. Kley & Lin 1992) and is easily distorted 
by even weak Lorentz forces. In this case we cannot ob- 
tain a closed equation because w* cannot be determined 
from a knowledge of the mass accretion rate. The essen- 
tial difficulty is that, while u r is the relevant mean velocity 
for mass accretion, u* is the appropriate mean for flux ac- 
cretion. These could be significantly disparate and could 
even differ in sign. The above discussion shows that the 
only solution to this problem is to solve explicitly for the 
vertical structure of the disk including the physics that 
determines the vertical profile of the radial velocity. 

4. MATHEMATICAL MODEL 

We now consider the set of equations governing the lo- 
cal vertical structure of a disk containing a mean poloidal 
magnetic field. The equations of Paper I will be augmented 
by the inclusion of additional physical effects. 



4.1. Basic equations for a thin disk 
The angular velocity is written as 

n = n {r) + n 1 (r,z,t), 

where 

1/2 



Q 



is the Keplerian value, and |Oi | -C O . 
The flux function is written as 

V> = i>o(r,t) + i>!(r,z,t), 

where |V>i| -C IV'ol, and we approximate 



(13) 



(14) 



(15) 



Z^TT> ^"Z-ST. (I 6 ) 



19Vo 



r dz 



r dr 



so that B z is independent of z. 

The required equations may be approximated as follows. 
The radial component of the equation of motion is 



— 2/9rfi fii 



B z dB r 
/j dz 



(17) 



where po is the permeability of free space. The azimuthal 
component is 



r dr 



B z dB 6 1 d 



p dz r 2 dr 



2 d£lQ 
dr 



(18) 



where v is the kinematic viscosity. The vertical component 
is 

= ^-1 ( P+ J| + ||), (19) 
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where p is the pressure. The poloidal part of the induction 
equation is 



chpo 
dt 



ru r B z 



( dB z dB r 



\ dr 



dz 



The toroidal part is 







dr 



B 



_9 

9x 



■5- »7 



<9z 



The energy equation is 



of ( dn y t; 

— = /ji/ r— — H 

az \ dr J Ho 



fdB r 



fdB 4 



where 



F = - 



\ dz 
16aT 3 dT 



\ dz 



(20) 



(21) 



(22) 



(23) 



3np dz 

is the radiative energy flux, with <r the Stefan-Boltzmann 
constant, T the temperature, and k the Rosseland mean 
opacity. 

These equations must be supplemented by constitutive 
relations specifying the equation of state, the opacity, the 
viscosity, and the magnetic diffusivity. We will adopt the 
ideal-gas equation of state, 



P 



kpT 
pm H ' 



(24) 



(where k is the Boltzmann constant, p the mean molecu- 
lar mass, and ma the mass of the hydrogen atom), and a 
generic power-law opacity, 



C K p x T y , 



(25) 



where C K , x, and y are constants. This includes the cases 
of Thomson scattering opacity (x = y = 0) and Kramers 
opacity (x = 1, y = —7/2). 

4.2. Turbulent viscosity and magnetic diffusivity 

The viscosity is less certain, and we will adopt the stan- 
dard prescription 

(26) 



ap 



where a is a dimensionless constant. For the magnetic 
diffusivity, we assume that the magnetic Prandtl number, 



Pm 



v 

j 

V 



(27) 



is constant. 

We will further assume either that a is a fixed parame- 
ter ('fixed alpha hypothesis'), or that a can adjust so as to 
keep the equilibrium at marginal magnetorotational stabil- 
ity, as explained in Section 4.8 below ('marginal stability 
hypothesis') . 

4.3. Neglected terms 

Several terms in the equations have been omitted on 
the grounds that the disk is thin and the solution should 
be stationary on the dynamical time-scale (although not 
necessarily on the viscous time-scale). The radial pressure 



gradient, the vertical variation of radial gravity, and the 
incrtial terms associated with the meridional flow, have 
been neglected as usual. However, enough terms have been 
retained to determine the profile of radial velocity in the 
absence of a magnetic field. 

Generally, fij has been neglected relative to Oo, and ipi 
relative to tpo, except where physically essential. It has 
also been assumed that the viscosity does not act on the 
shear components du r /dz and rdQ\/dz. Such terms are 
expected to be relatively unimportant in most cases of in- 
terest, and it was found that the inclusion of these terms 
increases the order of the differential system and makes it 
difficult to obtain a solution. It is especially unclear how 
to connect a viscous disk to an inviscid atmosphere, when 
these terms are included, without introducing an artificial 
boundary layer. 

We did not include a dynamo alpha-effect in the equa- 
tions for the mean magnetic field. Although there is no 
difficulty in principle in doing so, this effect is even less 
well understood than the turbulent viscosity and magnetic 
diffusivity, and we chose to avoid this additional compli- 
cation in the present study. 

Finally, in the energy equation ( p2[ ) it has been assumed 
that the heat generated in each annulus is radiated away 
locally and not advected through the disk. 

4.4. Local and non-local effects 

The equations of Section 4.1 describe the local vertical 
structure of an accretion disk with a mean poloidal mag- 
netic field. When supplemented with appropriate bound- 
ary conditions, as described below, they constitute a prob- 
lem similar in type to a stellar structure calculation, al- 
though the details are of course very different. Obviously it 
is possible in principle to include a more detailed equation 
of state, accurate opacity tables, a more sophisticated ap- 
proach to radiative transfer, and further possibilities such 
as convective energy transport. However, the principal un- 
certainties concern the viscosity and magnetic diffusivity. 

In fact, these equations are not strictly local in radius 
and time because some radial derivatives and one time- 
derivative remain. The case of dflo/dr = — 3flo/2r is triv- 
ial. The terms dipo/dt and dB z /dr in the poloidal part of 
the induction equation are retained in order to determine 
whether the magnetic flux migrates inwards (dipo/dt > 0) 
or outwards as a result of the local disk solution. The gra- 
dient dB z /dr can affect this result because it contributes 
to the diffusion of flux. Therefore dB z /dr appears as a pa- 
rameter of the local model and dipo/dt as an eigenvalue; 
note that both are independent of z. 

The case of d(pv)/dr is more problematic. This viscous 
term is retained in the angular momentum equation be- 
cause it partially determines the radial velocity, which in 
turn causes radial advection of magnetic flux and also af- 
fects the shape of the field lines, at least when the field is 
weak. This term may be approximated by arguing that 



pv * IT f (I) 



(28) 



in the neighborhood of the radius under consideration, 
where 

i?£ = / pv dz (29) 

J-H 



6 



Launching of jets 



is the vertically integrated dynamic viscosity, H the semi- 
thickness, and / an undetermined dimcnsionlcss function 
Under this assumption, 

dlnH 



dln(pu) _ dln(PE) 



d In? 



dlxi'i 



1 



d\n(pv) 

dhiz 



91m 



(30) 



The vertical derivative is available as part of the local so- 
lution, while 91n(PS)/91nr and d\nH/d\nr appear as 
additional dimensionless parameters, which can be esti- 
mated from the well-known behavior of the steady, non- 
magnetized solution (Shakura & Sunyaev 1973). 

In the limit of a non-magnetized disk, this prescription 
predicts the radial velocity in the disk as 



2r~ 



1 



<91nr 



1 



d In (pi/) 
dhiz 



dkiH 
dhir 



(31) 



4.5. Boundary conditions 



We have derived a sixth-order system of nonlinear ordi- 
nary differential equations (ODEs) in z. The dependent 
variables may be taken as p, T, F, fii, B r , and B^,. The 
variables p and u r are determined algebraically from these. 

The solution should be symmetric about the mid-plane, 
with 

F = B r = B^ = (32) 

at z = 0. 

The solution extends up to a photospheric surface z = H 
at which it is matched to a simple atmospheric model, dis- 
cussed below. The photospheric boundary conditions are 



F s = <jT 4 



and 



H 



npdz = — , 



(33) 



(34) 



where the subscript 's' denotes a surface value. We also 
have 

B r = B rs , £?0 = B<p s (35) 

there, where B rs and B^ are assigned values which, phys- 
ically, are determined by the solution of the global exterior 
outflow problem (not considered here). 

With these boundary conditions the equation of angular 
momentum conservation, obtained from a vertical integra- 
tion of equation (fish , has the form 



M d 



2irr 2 dr 



where 



M 



2B Z B<£ 
Po 



-lixr 



i_d_ 

r 2 dr 



vT,r 



3 d£lo 
dr 



H 



pu r dz 



(36) 



(37) 



H 



is the mass accretion rate, not necessarily constant. This 
shows the contributions to angular momentum transport 
from magnetic and viscous torques. 

The solution is determined as follows. One guesses the 
values of the quantities H, F s , and Oi B . This fixes the at- 
mospheric model and determines p s and T s . The quantities 
B rs and B^ are given. One must further guess dipo/dt to 
start the downward integration. The three guessed quan- 
tities should be adjusted to match the three symmetry 
conditions on z = 0. The main parameters of the model 
are then Q,q, E, B z , B rsi B$ s , a, and Pm. The quantity 
dtpo/dt is to be determined as an eigenvalue. 



4.6. Atmospheric model 

In the simplest atmospheric model, T, F, B ri and B^ are 
independent of z between the photosphere and the sonic 
point, with 



T = T S = T cB , 
F = F S = aT s \ 

B r = B rs = B z tan i, 

Ba, = B^, 



(38) 
(39) 
(40) 
(41) 



where i is the inclination of the field lines to the verti- 
cal. The constancy of F relies on there being little or no 
dissipation in the atmosphere, while the constancy of B 
relies on the atmosphere being magnetically dominated so 
that field line bending cannot be supported. That is, the 
plasma beta based on the poloidal magnetic field strength, 



2poP 
Bl + Bl : 



(42) 



should be less than unity. 

In the atmosphere, all diffusive terms in the equations 
are neglected. If dtpo/dt ^ 0, there must be a slow drift or- 
thogonal to the poloidal field lines, but otherwise the flow 
is constrained to follow the field. We know from Paper I 
that there will be a transonic outflow if i > 30° , otherwise 
a modified hydrostatic atmosphere. 

As in Paper I, the angular velocity in the atmosphere is 
determined by isorotation, 



3fi 

Qi = Qis + — — (z — H) tan i. 
2r 



(43) 



Again, this assumes that the atmosphere is magnetically 
dominated. The centrifugal-gravitational potential for 
matter constrained to follow the field is 



$ cg = -i(3tan 2 z - l)ng(z - z sonic ) 2 , (44) 



where 



= H 



1 + (3 tan 2 i - 1)" 



fti s r 
1 — 2 tan i 

n H 



(45) 



is the height of the sonic point in the case i > 30°. 

The density scale-height h s at the photosphere is given 

by 



nifll l-2^tani ) . (4<>) 



h s \ dz 

where c s is the isothermal sound speed, given by 



(47) 
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To a good approximation, the photospheric boundary con- 
dition (EmI) equates to 



exp 



-(1 



z-H 



This may be rearranged into the form 

1 + x\ / /iTOn \ / 16(7 



d(z-H). 
(48) 



ol +x T}+y = 



8fT 



V k J \3C K 
Q\ s r 



(49) 



Finally, in the case of a transonic outflow, the Mach 
number at the surface is determined from the equation 
(see Paper I) 



\(Ml 



1 )-\nM s = 



where 



(HqH — 20i s rtani) 



(50) 



(51) 



2(3 tan 2 i - 1) 

is the potential barrier to the outflow. The vertical mass 
flux density in the outflow is then 



m w = pu z — j\4 s p s c s cosi. 
4.7. Non-dimensionalization 



(52) 



We now rewrite the equations in a non-dimensional form 
suitable for numerical analysis. Given the local surface 
density E, the angular velocity Qq, and the constants ap- 
pearing in the constitutive relations, we identify 



jj h _ ^(2+x)/(e+x-2y)Q-{5-2y)/(6+x-2y) 



-(4-y)/(e+x-2y) 



16a 



-l/(6+x-2y) 



(53) 



as a natural unit for the semi-thickness of the disk. Nat- 
ural units for other physical quantities follow according 
to 

l2, 



U p = YXJ 



u T = nl (^M) u H 



U B 



U p — ~£Q Uh, 

u F = n u H u p , 

(poUp) 1 / 2 . 



(54) 
(55) 
(56) 



Note that the above expression for Uh can be obtained 
from the condition 



U F = 



16o-Z7| 



U r 



3C K U-U^U P U H 



(57) 



which is a dimensional analysis of the definition (|23|) of the 
radiative flux. 

There are two small dimensionless parameters in the 
problem, 

e = — 5= *i4- (58) 

Evidently e is a characteristic measure of the angular semi- 
thickness of the disk, while 6 is an inverse measure of the 
total optical thickness. 



Non-dimensional variables are then introduced using the 
transformations 



z = z* U H , 


H = 


H*U H , 


(59) 


p = p* u pi 


P = 


p* u p 


(60) 


T = T* U T , 


F = 


F*U F , 


(61) 


v = 14 Ct U H , 


V = 


i]* n u H , 


(62) 


u r = u r * ^IqUh, 






(63) 


B r = B r1f Ub, 


Brj, 


= Ub, 


(64) 



B z = B Z3f U B , ^- = ^m U H U B . (65) 
ot 



The transformed equations are 

— 2p*fii Jl , = B z *B' r% , 



(66) 



13 3 

= B Z *B'^ - -e(l + 2D u y,)p*v* + -jtU>n{P* v * z *) '> 

(67) 

= -p*z* -K - B r *B^ - BfaB'^, (68) 
V>* + u„B„ = 77* {eD B B z * - B' rt ) , (69) 

(70) 

F[ ='-!>. ^ -rih[( B r*) 2 + (Z'<f,*) 2 ], (71) 

T?- y 



o = --B„ + B w ai, + (i),y 
4' 



-T' 



(72) 
(73) 



where 



dln(vE) dlnff _ dlnB z 
£vs = — ?n 1 -Off = — , Df 



91nr 



<91nr 



<91nr 

are three dimensionless parameters, and the prime denotes 
differentiation with respect to z*. 

The dimensionless viscosity and magnetic diffusivity are 
given by 

v* = Pmr), = aT*. (75) 
The mid-plane symmetry conditions are 

F*(0) = 5 r *(0)=^*(0) = 0. (76) 

The photospheric boundary conditions are 

8F.(H.) = [T„{H*)} 4 (77) 

and 

[p*(H*)] 1+x [T*(H*)} 1+y = i*(l + x)[^-2n 1 *(^)tani], 

(78) 

together with 
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Finally, the definition of £ implies the normalization con- 
dition 

2 / p*dz» = 1. (80) 
Jo 

When the solution is found, the optical depth at the mid- 
plane can be determined from 



H 



Kpdz = 



1G 
35 



H, 



-X+x 



TUdz* 



(81) 



The fact that the small parameter e cannot be fully 
scaled out of the equations indicates that a strictly self- 
consistent thin-disk asymptotic solution cannot be ob- 
tained when all the physical effects we have considered 
are included. The e terms in equation d67j ) are expected 
to be important only in the case of a weak magnetic field. 
The e term in equation ( |69| ) may be important when the 
field is nearly vertical. 

4.8. Stability considerations 

Under the fixed alpha hypothesis the magnitudes of 
the turbulent viscosity and magnetic diffusivity are de- 
termined by the prescribed value of a. A problem with 
this approach, as will be seen in the next section, is that 
equilibrium solutions can then be obtained that are mag- 
netorotationally unstable even when the effect of the tur- 
bulent diffusivity on stability is taken into account. This 
suggests that the model may be physically inconsistent for 
these examples, because the 'channel solution' (Hawley & 
Balbus 1991) would continue to grow exponentially. Fur- 
thermore, there is a close connection between the stability 
of the equilibria and properties such as the shape of the 
field lines (Ogilvie 1998). We will find solutions that we 
believe to be unstable, in which the field lines bend more 
than once, and which have other undesirable properties. 

A resolution of these difficulties is suggested by the 
numerical simulations of magnetorotational turbulence, 
which indicate that the turbulence is much more vigorous 
in the presence of a mean poloidal magnetic field (Haw- 
ley, Gammie, & Balbus 1995), provided that the field is 
not so strong as to suppress the instability. Indeed, Stone 
et al. (1996) found it impossible to run a simulation in 
a stratified disk with a net vertical field, while Miller & 
Stone (2000) found a dramatic difference between simula- 
tions with no net field and those with a fairly weak uniform 
vertical field. 

In an effort to understand and model this behavior, we 
propose the following physical hypothesis for turbulent 
disks with a mean field: the value of a adjusts so that the 
equilibrium is marginally stable to the magnetorotational 
instability when the turbulent magnetic diffusivity is taken 
into account. We remark that a similar situation occurs in 
stellar convection, where an equivalent hypothesis can be 
used as a basis for the mixing-length theory. In a further 
example of this approach, Kippenhahn & Thomas (1978) 
modeled the outcome of a shear instability by imposing 
marginal stability according to the Richardson criterion. 

When the mean poloidal magnetic field is very weak, the 
instability favors small length scales and only a small value 
of a is required to suppress it. For stronger fields, the re- 
quired value of a will be larger. For even stronger fields 
with [3 ss 1 on the mid-plane, the instability is suppressed 



or nearly so even without any diffusivity, and a will be 
small or zero. This variation of a with the strength of the 
mean poloidal field is qualitatively in agreement with that 
found in numerical simulations (Hawley, Gammie, & Bal- 
bus 1995; Brandenburg 1998), suggesting that our physical 
hypothesis is reasonable. However, our model will require 
this hypothesis to extend into a parameter range that can- 
not be (or at least has not been) reproduced in numerical 
simulations of stratified disks. 

Although a magnetorotationally unstable disk may con- 
tain numerous unstable modes, it has been argued by 
Ogilvie (1998) that, in ideal MHD, the last mode to be 
stabilized is an axisymmctric mode with vanishing radial 
wavenumber (i.e. d/dr ~ r _1 rather than d/dr ~ H^ 1 ), 
and is the first mode of odd symmetry about the mid- 
plane. Marginal stability of the annulus is then imposed 
by locating the marginal stability condition for this crit- 
ical mode. We conjecture that this remains true when 
dissipation is included. 

From a consideration of the linearized equations in the 
limit of vanishing eigenfrequency, we find that such a mode 
should involve only horizontal motions and satisfy the hor- 
izontal components of the perturbed equation of motion, 



i/of^o Au r = 
2 p 



B z d AB r 
/io dz 

B z d ABa, 



dz 



and the perturbed induction equation, 







= B z 



AB r 



d Au r 



d 

dz {1] - 



d AB, 



dz 



B 



dAu 4 
dz 



d ( d AB d 



dz^- 



dz 



(82) 
(83) 

(84) 
(85) 



where A denotes a linearized Eulerian perturbation. These 
equations should be compared with the unperturbed equa- 
tions (|l|)-(|l|) and Note that this problem is 
different from the marginal stability problem considered 
in Paper I, where the magnetic diffusivity was zero. 

The critical mode is expected to be the first mode of 
odd symmetry about the mid-plane (Ogilvie 1998). This 
should satisfy the symmetry conditions Au r — Au^ = 
on the mid-plane, and the photospheric boundary condi- 
tions AB r = AB,f, = (since B rs and B^ s are fixed in the 
local model). For such a mode, the linearized equations 
may be combined into the second-order ODE 



d_ 

dz 



[Bl 



fi PV 2 ^l\ d AB 



\p p B 2 Z 
or, in dimensionless form, 



dz 



3f^AB r = 0, (86) 



P* 



^ + ?^k\ AD', 



3 AB„ = 0, 



(87) 



with boundary conditions 

AB;„(0) = AB„(H*) = 0. 



We are now in a position to determine the value of a 
as follows. For any given value of a and suitable values of 
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the other parameters, we might expect to find an equilib- 
rium solution. However, the equilibrium will not in general 
possess a marginal mode satisfying the correct boundary 
conditions. By integrating the equations for a marginal 
mode simultaneously with those defining the equilibrium, 
we aim to tune a until a marginal mode is obtained. If 
this mode is the first mode of odd symmetry, we then be- 
lieve we have determined the value of a corresponding to 
a marginally stable equilibrium. 

A crude estimate of this may be obtained by approxi- 
mating the second derivative with respect to as —H~ 2 . 
This leads to the estimate 

Vl « ^ ( ZHl - ^) , (89) 
P* \ P* J 

which has the properties described above. 

It might be argued that our marginal stability condi- 
tion is too strong because the presence of one or two un- 
stable modes of relatively long radial wavelength might 
not contribute significantly to turbulent transport. How- 
ever, if the critical mode is not stabilized, it will lead to 
the growth of the fundamental 'channel solution', which 
is very efficient in transporting angular momentum and 
probably hydrodynamically unstable in three dimensions, 
leading to enhanced turbulence. Admittedly, the descrip- 
tion of the disk under these circumstances is a matter of 
some uncertainty. 

5. NUMERICAL INVESTIGATION 

5.1. Numerical method 

The dimensionless ODEs are integrated from the pho- 
tosphere Z.J, = to the mid-plane = 0. The depen- 
dent variables are p*, T*, F*, ill*, B r „, and B^. The 

values of H*, fii*(if*), and ip* are guessed and 

then adjusted by Newton-Raphson iteration to match the 
symmetry conditions on the mid-plane. In the marginal 
stability model, the equations for the marginal mode are 
integrated simultaneously. An arbitrary normalization 
AB' r3f (H*) = 1 of the linear problem is applied, and 
the value of a is also guessed and adjusted by Newton- 
Raphson iteration to match the symmetry condition for 
the mode on the mid-plane. 

To obtain the derivatives B'^ and fi^ some algebra is 
required. Eliminating it r „, p'^, and B' r3t we obtain 



* 2PmB^(2B zr - ZeDnaB^z*) 
x | - 2Pm B z *ip* + 4ap*Sli* 

+e [2D B + 3Pm (1 + 2D vYl - 2D H )] aB^T, 
+6PmeD H aB„ (B z *zl - 2B r *n u z*) j. (90) 

By multiplying this expression by 77* and differentiating, 
we find Qj^ from equation (jTOj) . 

5.2. Unmagnetized solution 

In the absence of a mean magnetic field, a solution is 
obtained by omitting the induction equation and setting 



B = elsewhere. This also implies Q\ = 0. The form 
of the solution for p*, T*, and depends only on the 
dimensionless parameter 8 (for a given opacity law), al- 
though the variables also have power-law scalings with a. 
If u r * is required, there are further dependences on D v y. 
and Dh, and the value of u r * is proportional to e as a con- 
sequence of the scalings we have adopted for the general 
problem. 

We focus on the case of Thomson scattering opacity 
(x — y — 0, C K w 0.33cm 2 g" 1 ). In a steady disk, 
far from the inner edge, we may then set D v y, = and 
D H = 21/20 (Shakura & Sunyaev 1973). We also take 
a = 0.1. For the purposes of illustration, we consider 
a location at 1000 Schwarzschild radii from a black hole 
of mass IOM0, i.e. r = 2.95 x 10 9 cm. For a surface 
density £ = 10 gem , and assuming fi — 0.6, we find 
illustrative values Uh = 6.34 x 10 7 cm, U p — 1.58 x 
10- 4 gcm- 3 , U p = 3.27 x 10 10 dyn, U T = 1-50 x 10 6 K, 
U F = 4.70 x 10 17 ergcm- 2 s-\ and U B = 6.41 x 10 5 G. 
Then e = 2.15 x 10~ 2 and 8 = 1.62 x 10~ 3 . 

The numerically determined unmagnetized solution is 
shown in Fig. 1. It has a dimensionless photospheric height 
of — 1.70 and an optical depth at the mid-plane of 
r c = I + §5 _1 . Thus our illustrative values correspond to 
H/r = 0.0364 and r c = 1650, representing a geometrically 
thin and optically thick disk. The illustrative accretion 
rate is M = 3.55 x 10 18 gs~ 1 , which is less than the Ed- 
dington accretion rate 1.69 x 10 19 gs _1 for an accretion 
efficiency of 0.1. 

The profile of radial velocity is of particular interest, 
since this will affect the advection of the magnetic field. 
The radial velocity is positive on the mid-plane and be- 
comes negative at larger z, the density- weighted average 
being of course negative. This result has been found pre- 
viously by several authors (e.g. Urpin 1984; Kley & Lin 
1992, where further analysis of this phenomenon may be 
found). In the present model it results from the vertical 
profile of the viscosity, specifically from the dlnH/dlnr 
term in equation (|3l]) . 

5.3. Magnetized solutions under the fixed alpha 
hypothesis 

We first adopt the fixed alpha hypothesis and search for 
magnetized solutions, taking parameter values a = 0.1, 
-B^s* = 0, Pm = 1, D B = in addition to x = y = 0, 
D„£ = 0, Dh = 21/20, and the illustrative values for e 
and 8. There are two limits in which physically acceptable 
solutions are easily obtained over a wide range of incli- 
nation angles: very weak fields, and strong fields. We 
illustrate two such cases, with i = 45°, in Figs 2 and 3. 

In the very weakly magnetized solutions, the magnetic 
field behaves kincmatically and the equilibrium structure 
of the unmagnetized solution is not significantly distorted 
(e.g. the radial velocity profile is almost indistinguishable 
from that in the unmagnetized solution). The poloidal 
field lines bend smoothly according to a balance between 
advection and diffusion. Isorotation is not enforced (there 
is only a minuscule deviation from Keplerian rotation) and 
the field is strongly wound up (\B^,\ > \B Z \) in the exam- 
ple shown, with i = 45°. The poloidal field lines continue 
to bend as the photosphere is approached. This is permit- 
ted because f3 is very large throughout. Our atmospheric 
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model, which assumed (3 < 1 , is not relevant here and cer- 
tainly the outflow solution should be disregarded. It is not 
clear how to treat the atmosphere in a case such a this, 
because it is not understood whether the radial flow and 
magnetic diffusivity continue above the photosphere. 

By setting _B^ S = 0, we have ensured that the radial 
mass flux is caused by viscosity alone (cf. eq. |36|). As 
expected, the radial velocity causes the flux to migrate 
inwards for almost vertical fields, but outwards for larger 
inclination angles. With Pm = 1 and Db = 0, one finds 
the intermediate case ip* — at i — 0.485 = 27.8°. This 
is much larger than H/r = 0.0364, indicating that our 
concerns expressed in Section 3 were well founded. Alter- 
natively, one can prevent the flux from being expelled (i.e. 
achieve ^* = 0) when i = 45° by increasing Pm from 1 to 
1.90. Increasing Db does not change ip* significantly. 

For strong fields, we recover solutions very similar to 
those we obtained in Paper I. In the example shown in 
Fig. 3, the disk is significantly compressed by the Lorentz 
force (compare the plot of p„ with the unmagnetized case) , 
and is hotter (cf. T*) and more luminous (cf. F*) than its 
unmagnetized counterpart. Most of the bending of the 
poloidal field lines occurs near the mid-plane where (3 is 
close to unity. The field enforces isorotation, resulting in 
sub-Keplerian rotation, and very little toroidal field is pro- 
duced. 

A feature of the solution is that the radial velocity, 
although everywhere subsonic, is rather large and non- 
uniform in direction. The reason for the profile of IS 
that the flux is being expelled rapidly (ift* — —0.551) by 
turbulent diffusion. In the upper layers, where (3 is small, 
the fluid must follow the field and therefore flows rapidly 
outwards. To achieve the net accretion rate imposed by 
the vertically integrated viscous stress, a strong inflow is 
then required close to the mid-plane. Such a profile might 
even be considered advantageous for jet launching, since 
the outflow itself must involve a positive radial velocity. 

With B z * — 2, as in this example, and assuming Pm = 
1, zero flux migration is found at i = 0.0346 = 1.98°. This 
is now close to H/r = 0.0364 as expected from simple ar- 
guments, probably because the field-line bending occurs 
close to the mid-plane. 

On the other hand, these considerations may be incon- 
sistent because the example in Fig. 3 is magnetorotation- 
ally stable and therefore unlikely to be turbulent. Of 
course, if there is no turbulent magnetic diffusivity, the 
origin of the viscosity should also be questioned. 

For a range of intermediate field strengths, the situation 
is more complicated. Solutions are not found in which the 
field lines bend in the simple way seen in Figs 2 and 3. 
Instead, branches of solutions appear in which the field 
lines bend several times as they pass through the disk. An 
example is shown in Fig. 4. 

According to Ogilvie (1998, Theorem 1) such a solution 
would be magnetorotationally unstable in ideal MHD. The 
multiple bending is closely related to the appearance of the 
'channel solution' which is the first stage of the magnetoro- 
tational instability for a vertical field. We conjecture that 
such multiple-bending solutions are also unstable when 
the turbulent magnetic diffusivity is taken into account 
(in both the equilibrium and perturbation equations), al- 
though this has not been proven. If correct, this would im- 



ply that solutions of this type are physically inconsistent 
and ought to be disregarded, which in turn would suggest 
that no steady solution is possible for these intermediate 
field strengths. 

However, our alternative proposal, the marginal stability 
hypothesis, allows for a possible solution of this difficulty. 

5.4. Magnetized solutions under the marginal stability 
hypothesis 

Under the marginal stability hypothesis, the strength of 
the turbulence (quantified through the parameter a) is just 
sufficient to bring the equilibrium to marginal magnetoro- 
tational stability. This enables us to find single-bending 
solutions in a continuous range of field strengths from very 
small values up to the strength at which the instability is 
suppressed even without a turbulent diffusivity. We adopt 
the same parameters as in the previous section, except 
that a is now to be determined as part of the solution. 
Fig. 5 shows how the required value of a varies with B z * 
for solutions with vertical fields (i = 0). 

This has the general form expected from the crude es- 
timate, equation (|89|). However, it is somewhat worrying 
that values of a in excess of unity may be required. It is 
often argued that a should not exceed unity because the 
turbulence would then have to be supersonic, or to have 
a correlation length greater than the disk thickness (e.g. 
Pringle 1981). It is not clear whether these constraints 
necessarily apply to magnetorotational turbulence, which 
is dominated by anisotropic Maxwell stresses whose cor- 
relation length in the azimuthal direction could exceed H 
(e.g. Armitage 1998), but which could be limited instead 
by magnetic buoyancy. Moreover, the effective transport 
coefficients may themselves be anisotropic in the presence 
of a strong mean field (e.g. Matthaeus et al. 1998) and the 
magnetic Prandtl number may also differ from unity. In- 
deed, it is not certain whether the effects of the turbulence 
can be described adequately in terms of an effective dif- 
fusivity; it may be that coherent magnetic structures are 
formed. Nevertheless, the question remains as to whether 
the turbulence can truly reach a level sufficient to achieve 
marginal stability in the presence of a significant mean 
field. The numerical simulations by Hawley et al. (1995) 
are in qualitative agreement with Fig. 5 but suggest that, 
when the field strength is just below the stability bound- 
ary for a given computational domain, the channel solu- 
tion may continue to grow without degenerating into tur- 
bulence. However, simulations of stratified models with 
fairly strong mean vertical fields do not appear to have 
been successful. 

An example solution calculated under the marginal sta- 
bility hypothesis is shown in Fig. 6. This has a — 3.81 and 
should be contrasted with the solution in Fig. 4 which has 
the same vertical field strength but a fixed a = 0.1. The 
marginally stable solution has field lines with a single bend 
which become straight not far below the photosphere. The 
pressure and temperature decline smoothly and monoton- 
ically with increasing height. The field is not significantly 
wound up (l-B^I < \B Z \). 

5.5. Mass loss rates for jet-launching disks 

The principal aim of this paper was to determine how 
the mass loss rate in the outflow, m w , varies with the 
strength and inclination of the magnetic field. In Fig. 7 



Ogilvie & Livio 



11 



we show the result of this calculation under the fixed alpha 
hypothesis with a — 0.1 and other parameters as given in 
Section 5.3. Note that the reciprocal of m w /(£fi) is ap- 
proximately the number of orbits in which the disk would 
be evaporated if not replenished by the accretion flow. In 
Fig. 8 we plot the quantity 47rr 2 m w /M as a function of 
B„ for solutions with i = 45°. This alternative dimen- 
sionless measure of the outflow rate is the local mass loss 
rate per unit logarithmic interval in radius, divided by the 
accretion rate. 

The outcome is as we found in Paper I. The solutions 
shown arc magnetorotationally stable. Close to the edge 
of the solution manifold, the outflow is very vigorous, but 
if the field strength is increased by only a factor of two 
above the stability boundary, the outflow is suppressed by 
twenty orders of magnitude or so. As explained in Paper I, 
this happens because the disk becomes significantly sub- 
Keplerian (cf. Fig. 3) and the outflow experiences a large 
potential barrier (whose strength is roughly proportional 
to Bf). This means that, in the absence of additional 
heating, external irradiation, or other driving mechanisms, 
outflows are suppressed from strongly magnetized disks. 
For a fixed field strength, the outflow is maximized at an 
intermediate inclination angle of 40-50°. 

The corresponding results for the marginal stability hy- 
pothesis are shown in Figs 9 and 10. Note that this is a 
complementary region of parameter space corresponding 
to turbulent disks. The behavior of the mass loss rate is 
now completely different and perhaps more intuitive: it 
increases monotonically with increasing B zst and with in- 
creasing i. However, the solution branch terminates before 
excessive mass loss rates are achieved. This suggests that 
a large potential barrier is not incurred for turbulent disks 
and that they are more promising for jet launching. The 
potential barrier is smaller because the magnetic field is 
weaker, resulting in a smaller Lorentz force and a smaller 
deviation from Kcplerian rotation. 

These two sets of results were obtained under differ- 
ent physical assumptions. Nevertheless, they both suggest 
that efficient jet-launching solutions are found in a limited 
range of field strengths, and in a limited range of inclina- 
tion angles in excess of 30°. In both cases there are dif- 
ficulties in interpreting the solutions. The more strongly 
magnetized solutions obtained under the fixed alpha hy- 
pothesis are magnetorotationally stable, and the origin of 
the dissipation in the disk remains unclear. The more 
weakly magnetized solutions obtained under the marginal 
stability hypothesis typically require values of a in excess 
of unity to bring the equilibrium to marginal magnetoro- 
tational stability. 

6. DISCUSSION 

We have developed a model of the local vertical struc- 
ture of magnetized accretion disks that launch magneto- 
centrifugal outflows. Given certain assumptions concern- 
ing the dissipative processes (turbulent viscosity and mag- 
netic diffusivity) in the disk, we have shown that it is pos- 
sible to compute the mass loss rate in the outflow as a 
function of the surface density and the strength and incli- 
nation of the poloidal magnetic field. The net accretion 
rates of mass and magnetic flux are also determined. This 
information is precisely complementary to that obtained 



from numerical simulations of the subsequent acceleration 
and collimation of a 'cold' outflow (e.g. Krasnopolsky et al. 
1999). 

The following result appears to be quite robust. For 
disks in which the mean poloidal magnetic field is suffi- 
ciently strong to stabilize the equilibrium against the mag- 
netorotational instability, we find that the mass loss rate 
decreases extremely rapidly with increasing field strength, 
and is maximized at an inclination angle of 40-50°. For 
turbulent disks with weaker mean fields, we find that 
the mass loss rate increases monotonically with increas- 
ing strength and inclination of the field, but the solu- 
tion branch terminates before excessive mass loss rates are 
achieved. This suggests that turbulent disks with moder- 
ate mean fields are more promising for jet launching, but 
there may be situations in which a steady solution is im- 
possible. 

For each solution we have determined the net rate of 
flux migration, which depends on a competition between 
inward dragging by the accretion flow, and outward trans- 
port due to turbulent diffusion. The results depend on 
the effective magnetic Prandtl number of the turbulence, 
which has never been measured. However, we find that in- 
ward migration is more likely to occur in the case of a weak 
mean field, which bends at greater heights in the disk. In 
this case, inward dragging may be many times more ef- 
fective than estimated in previous studies (Lubow et al. 
1994a), with the result that a net inward migration of flux 
may occur even when the inclination angle is sufficient for 
jet launching, provided that the disk is not very thin. This 
issue requires further analysis and a better understand- 
ing of the kinematic behavior of the magnetic field. For 
stronger fields, we find fairly good agreement with previ- 
ous studies, and the jet-launching configurations are much 
more difficult to maintain against dissipation. We specu- 
late that this may lead naturally to a situation in which 
the flux is regulated to a value suitable for jet launching. 
However, we have not modeled the possible effect of an 
internal dynamo in the disk. In addition, it is quite prob- 
able that instabilities will arise when the global coupling 
between the outflow and the flux evolution is considered 
(Lubow, Papaloizou, & Pringlc 1994b; Cao & Spruit, in 
preparation) . 

Several other authors have considered the vertical struc- 
ture of magnetized disks and the problem of the disk-jet 
connection. Our formulation is distinctive in that we have 
shown how to calculate the mass loss rate and the accre- 
tion rates of mass and magnetic flux at any radius from a 
knowledge of the local surface density and the strength and 
inclination of the magnetic field. When coupled with a nu- 
merical simulation of the outflow beyond the sonic point, 
this would provide a closed evolutionary scheme for a time- 
dependent magnetized accretion disk with an outflow. It 
is crucial that we allowed the external magnetic torque on 
the disk to be determined by the exterior outflow solu- 
tion rather than the interior solution of the vertical disk 
structure. 

There are further differences in the detail of our ap- 
proach. Wardle & Konigl (1993) discussed many of the 
issues that we have addressed, but we differ from them in 
considering an optically thick disk with turbulent trans- 
port and energy dissipation, as opposed to an isothermal, 
inviscid disk with ambipolar diffusion. We are also able to 
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calculate the net rate of flux accretion rather than specify- 
ing it as a parameter. Similar comparisons may be made 
with the analysis of Li (1995). The model of Campbell 
(1999) shares some features with the present paper, but his 
outflow solution treats the magnetic field lines as parabo- 
lae throughout the disc and the transonic region, while 
our solutions (e.g. Fig. 6) indicate that the field lines be- 
come straight in the atmosphere, provided that the plasma 
beta is less than unity. If (3 > 1 in the atmosphere, the 
magnetocentrifugal mechanism does not work. 

Interesting comparisons may be made with the steady, 
self-similar model of Casse & Ferreira (2000). In common 
with those authors, we find that values of a exceeding 
unity, and/or anisotropic effective turbulent diffusion coef- 
ficients, may be helpful or necessary in obtaining plausible 
solutions for efficient jet-launching disks. 

Admittedly, our model has some limitations. Our sim- 
ple modeling of the effects of the turbulence is consistent 
with the very limited information available from numeri- 
cal simulations, although we have not included a dynamo 
a-effect. We have assumed that the atmosphere is magnet- 
ically dominated so that the poloidal field lines act as rigid 
channels for the outflow. This constraint could be relaxed 
in future work at the expense of introducing considerable 
complications in solving for the atmospheric flow, and in 



matching to the exterior solution. 

Future work should aim at obtaining a better under- 
standing of the interaction of the mean magnetic field with 
the turbulence. Numerical simulations might be used to 
measure the effective magnetic diffusivity tensor and how 
it depends on the strength of the field. A crucial issue 
is whether, when the mean field is a significant fraction 
of the value required for magnetorotational stability, the 
'channel solution' of the instability persists or degenerates 
into turbulence. 

Ideally the whole problem that we have defined, includ- 
ing the correct boundary conditions, would be solved using 
a numerical simulation. We have included some non-local 
effects (Section 4.4) but most of the relevant terms are in- 
cluded in the shearing-box local accretion disk model (e.g. 
Stone ct al. 1996). It may be difficult, however, to resolve 
the disk adequately from the mid-plane to the sonic point, 
although the techniques used by Miller & Stone (2000) 
should be helpful in achieving this. 

GIO acknowledges support from the STScI visitor pro- 
gram and from Clare College, Cambridge. ML acknowl- 
edges support from NASA Grants NAG5-6857 and GO- 
07378. We thank Jim Pringle and Henk Spruit for helpful 
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Fig. 1. — Profiles of pressure, temperature, radiative energy flux, and radial velocity for an unmagnetized model with Thomson opacity. 
Note that in this case the radial velocity is directed outwards on the mid-plane and becomes inward at greater height. 
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Fig. 2. — Profiles of pressure, temperature, radiative energy flux, radial velocity, deviation from Kcplerian angular velocity, radial magnetic 
field, toroidal magnetic field, and poloidal plasma beta, and shape of the poloidal field lines, for a very weakly magnetized model with 
B zt = 10~ 4 (i.e. B z = 64.1 G for our illustrative parameters) and i = 45°. 
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Fig. 3. — A strongly magnetized model with B z * = 2 (i.e. B z = 1.28 X 10 6 G for our illustrative parameters) and i = 45°. 
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Fig. 4. — A peculiar model with B zf =0.1 (i.e. B z = 6.41 X 10 4 G for our illustrative parameters) and i = 30°. Note the multiple bending 
of the poloidal field lines, which probably indicates a magnctorotationally unstable configuration, and other undesirable characteristics. 
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Fig. 5. — Variation of the viscosity parameter a with the mean vertical field strength B zt under the marginal stability hypothesis, for 
solutions with vertical fields (i = 0). 
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Fig. 7. — Mass loss rate as a function of the strength and inclination of the magnetic field, under the fixed alpha hypothesis. Contours of 
log 10 (m w /ESl) arc plotted, from —14 (darkest) to —4 (lightest) with unit spacing. The solutions continue to larger values of B zt , but do not 
extend below the solid line. 




Fig. 8. — Ratio of the mass loss rate to the accretion rate as a function of the strength of the magnetic field for equilibria with i = 45°, 
under the fixed alpha hypothesis. The same solutions appear in Fig. 7. 
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Fig. 9. — Mass loss rate under the marginal stability hypothesis. Contours of log 1 Q(rn w /Sf2) are plotted, from —14 (darkest) to —4 
(lightest). The solutions continue to smaller values of i, but do not extend to the right of the solid line. Below the dashed line, the atmosphere 
is not magnetically dominated (j3 > 1 at the photosphere) and the outflow solution should not be trusted. 
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Fig. 10. — Ratio of the mass loss rate to the accretion rate as a function of the strength of the magnetic field for equilibria with i = 35°, 
under the marginal stability hypothesis. The same solutions appear in Fig. 9. 



